/*
 * \file CplexStabSolver.cpp
 * \brief Implementación de la clase CplexStabSolver
 *
 *  Created on: Jan 24, 2014
 *      Author: ik
 */

#include <CplexStabSolver.h>
#include <Grafo.h>
#include <Desigualdad.h>
#include <SeparatorParameters.h>
#include <RDCliquer.h>
#include <HeuristicasGolosas.h>
#include <TimeLimit.h>
#include <cmath>

class RoundingHeuristicI;
class BranchingI;

//macro que se necesita para portabilidad, parece.
ILOSTLBEGIN

//macro para implementar fácilmente el callback para generar cortes.
ILOUSERCUTCALLBACK2(CutGenerator, CplexSeparator *, _separator, IloNumVarArray, _x)
{
	if (TimeLimit::getInstance().limitReached())
		return;

	CplexSolucion solucion(this, _x);

	list<Desigualdad *> desigualdades;
	int nnodes = getNnodes();
	int incumb = getIncumbentObjValue();
	int obj = getObjValue();
	if (nnodes == 0 || obj > 1.10*incumb)
	{
		if (nnodes == 0 || (nnodes % 100 == 0 && obj >= 1.25*incumb) || (nnodes % 50 == 0 && obj >= 1.35*incumb) || (nnodes % 25 == 0 && obj >= 1.50*incumb) || obj >= 1.75*incumb)
			desigualdades = _separator->separar(solucion);
		else {
			//	getEnv().out() << "Solucion de valor " << getObjValue() << " y best bound " << getBestObjValue() << "..." << endl;
			//	getEnv().out() << "La solución es [";
			//	for (int i = 0; i < _x.getSize(); i++)
			//		getEnv().out() << " (" << i << ")" << solucion[i];
			//	getEnv().out() << " ]" << endl;

			desigualdades = _separator->separarReserva(solucion);
		}
	}
	else
		_separator->getReservaCliques()->haceInactiva();

//	int ndesig = desigualdades.size();
//
//	if (ndesig == 0)
//		getEnv().out() << "Solucion de valor " << getObjValue() << " y best bound " << getBestObjValue() << "..." << endl;

	int nclique = 0;
	int nrango = 0;
	int ncorte = 0;
	while (!desigualdades.empty())
	{
		CplexDesigualdad * d = (CplexDesigualdad *) desigualdades.back();
		ncorte++;
		desigualdades.pop_back();
		if (d->isRank())
		{
			if (d->getRHS() == 1)
			{
				nclique++;
				if (nnodes > 0)
					add(*d, IloCplex::UseCutFilter);
				else
					add(*d);
			}
			else
			{
				nrango++;
				add(*d, IloCplex::UseCutFilter);
			}
		}
		else
			add(*d, IloCplex::UseCutFilter);
//			getEnv().out() << "-----> Corte con " << d->getSoporte().size() << " vertices: ";
//		if (nnodes > 0)
//			d->print();
		delete d;
	}
	//		if (ncorte > 0)
	//			getEnv().out() << "Hay " << ndesig << " desigualdades, incluyendo " << nrango << " de rango y " << nclique << " de clique" << endl;
}

//macro para implementar fácilmente el callback para generar cortes.
ILOHEURISTICCALLBACK3(RoundingHeuristic, HeuristicasGolosas *, _heuristica, IloNumVarArray, _x, Lista&, _intera)
{
	if (getNnodes() % 2 != 0) // || getNiterations() % 2 != 0)
		return;

	CplexSolucion solucion(this, _x);

	//	getEnv().out() << "Solucion de valor " << getObjValue() << " y incumbent " << getIncumbentObjValue() << "..." << endl;

	if (!TimeLimit::getInstance().limitReached())
	{
		_intera.clear();
		int c = _heuristica->rounding(solucion, _intera);
		if (c > getIncumbentObjValue()) {
			IloNumArray h(getEnv(), _x.getSize());
			for (int i = 0; i < _x.getSize(); i++)
				h[i] = 0;
			for (int i = 0; i < _intera.cantidad(); i++)
				h[_intera[i]] = 1;
			setSolution(_x, h, c);
		}
	}
}

class Fijos : public IloCplex::MIPCallbackI::NodeData
{
	VertexSet * _vertices;

	Fijos() : _vertices(NULL) {};
	Fijos(VertexSet * set) : _vertices(new VertexSet(*set)) {};
public:
	Fijos(int n) : _vertices(new VertexSet(n)) {};
	~Fijos() { if (_vertices != NULL) delete _vertices; };
	VertexSet& getVertices() { return *_vertices; };
	Fijos * newIstance() { Fijos * ret = new Fijos(); ret->_vertices = _vertices; _vertices = NULL;  return ret; };
	Fijos * newIstance(int v) { Fijos * ret = new Fijos(_vertices); ret->_vertices->add(v); return ret; };
};

//macro para implementar fácilmente el callback para ramificación.
ILOBRANCHCALLBACK2(Branching, CplexBranching *, branching, IloNumVarArray, _x)
{
	int n = _x.getSize();
	double max = 0.0;
	double min = 1.0;
	double x[n];
	int ix[n];
	int nx = 0;
	for (int i = 0; i < n; ++i) {
		double m = getValue(_x[i]);
		if (m < 1.0-1.0e-05 && m > 1.0e-05) {
			x[nx] = m;
			ix[nx++] = i;
			if (m > max)
				max = m;
			if (m < min)
				min = m;
		}
	}

	max = (max+min)/2.0;
//	max = max/2.0;
//	max = 0.5;
	int v = 0;
	double xv = 0.5;
	for (int i = 0; i < nx; ++i) {
		double nxv = fabs(x[i]-max);
		if (nxv < xv) {
			v = ix[i];
			xv = nxv;
		}
	}

	Fijos * fijos = (Fijos *) getNodeData();
	bool newFijos = (fijos == NULL);
	if (newFijos)
		fijos = new Fijos(n);

	list<Desigualdad *> desigualdades = branching->branchZero(v, fijos->getVertices());
	IloRangeArray constraints(getEnv());
	while (!desigualdades.empty())
	{
		CplexDesigualdad * d = (CplexDesigualdad *) desigualdades.back();
		desigualdades.pop_back();
		constraints.add(*d);
//		d->print();
		delete d;
	}
	constraints.add(_x[v] <= 0);
	makeBranch(constraints, getObjValue(), fijos->newIstance(v));

	desigualdades = branching->branchOne(v);
	if (!desigualdades.empty())
	{
		constraints.clear();
		while (!desigualdades.empty())
		{
			CplexDesigualdad * d = (CplexDesigualdad *) desigualdades.back();
			desigualdades.pop_back();
			constraints.add(*d);
			//		d->print();
			delete d;
		}
		constraints.add(_x[v] >= 1);
		makeBranch(constraints, getObjValue(), fijos->newIstance());
	}
	else
		makeBranch(_x[v], 1.0, IloCplex::BranchUp, getObjValue(), fijos->newIstance());

	if (newFijos)
		delete fijos;
}

CplexSolucion::CplexSolucion(IloCplex::ControlCallbackI * model, IloNumVarArray& x) : Solucion(x.getSize()), _model(model), _x(x) {
	_nfrac = 0;
	_nzeroes = 0;
	for (int i = 0; i < _n; i++) {
		if (operator [](i) < 1.0-1.0e-5) {
			_fracVertices.add(i);
			if (operator [](i) < 1.0e-5)
				_nzeroes++;
			else
				_nfrac++;
		}
	}
}

CplexSolucion::CplexSolucion(IloCplex::ControlCallbackI * model, IloNumVarArray& x, int n) : Solucion(n), _model(model), _x(x) {
}

inline const double CplexSolucion::operator [](int i) const
{
	return _model->getValue(_x[i]);
}

inline double CplexSolucion::operator [](int i)
{
	return _model->getValue(_x[i]);
}

inline void CplexSolucion::setValue(const int i, const double valor)
{
	throw new IloNotImplemented();
}

inline double CplexSolucion::value()
{
	return _model->getObjValue();
}

inline void CplexSolucion::print() const
{
	cout << "x = [";
	DoubleValueArray::print();
	cout << " ]" << endl;
}

DoubleValueArray* CplexStabSolver::solve(SeparatorParameters& parameters)
{
	DoubleValueArray * ret = 0;
	IloEnv env;
	try
	{
		IloModel model(env);

		int n = _grafo.getN();

//		for (VertexSetIterator it = _grafo.vertices().begin(); it != _grafo.vertices().end(); ++it)
//		{
//			if (_grafo.isStable(_grafo.vecinos(*it)))
//				cout << "vecindad de " << *it << " es independiente con " << _grafo.vecinos(*it).size() << " vertices!" << endl;
//			else if (_grafo.isInducedTree(_grafo.vecinos(*it)))
//				cout << "vecindad de " << *it << " es un árbol con " << _grafo.vecinos(*it).size() << " vertices!" << endl;
//		}
//		exit(0);

		IloEnv env = model.getEnv();

		//variable definition
		IloNumVarArray & x = *(new IloNumVarArray(env));
		for (int i = 0; i < n; ++i) {
			string name("x"+to_string(i));
			x.add(IloIntVar(env, 0, 1, name.c_str()));
		}

		//objective function
		IloNumExpr fobj(env);
		for (int i = 0; i < n; ++i)
			fobj += x[i];
		model.add(IloMaximize(env, fobj));

		IloCplex cplex(model);

		//callbacks
		CplexSeparator * separator = new CplexSeparator(_grafo, x, parameters);
		CplexBranching * branching = new CplexBranching(_grafo, x, separator->getReservaCliques());
		HeuristicasGolosas * heuristicas = new HeuristicasGolosas(_grafo);
		Lista * intera = new Lista(_grafo.getN());

		cplex.use(CutGenerator(env, separator, x));
		cplex.use(RoundingHeuristic(env, heuristicas, x, *intera));
		cplex.use(Branching(env, branching, x));

		//constraints
		IloRangeArray constraints(env);
		list<Desigualdad *> desigualdades = branching->branchNone();
		while (!desigualdades.empty())
		{
			CplexDesigualdad * d = (CplexDesigualdad *) desigualdades.back();
//			d->print();
			desigualdades.pop_back();
			constraints.add(*d);
			delete d;
		}
		model.add(constraints);

		//seteamos los parámetros
		cplex.setParam(IloCplex::NodeLim, 0);

		cplex.setParam(IloCplex::Cliques, -1);
		cplex.setParam(IloCplex::ZeroHalfCuts, -1);
		cplex.setParam(IloCplex::FracCuts, -1);
		cplex.setParam(IloCplex::ImplBd, -1);
		cplex.setParam(IloCplex::Covers, -1);
		cplex.setParam(IloCplex::DisjCuts, -1);
		cplex.setParam(IloCplex::FlowCovers, -1);
		cplex.setParam(IloCplex::FlowPaths, -1);
		cplex.setParam(IloCplex::GUBCovers, -1);
		cplex.setParam(IloCplex::MCFCuts, -1);
		cplex.setParam(IloCplex::MIRCuts, -1);
//		cplex.setParam(IloCplex::TiLim,TimeLimit::MAX_SECONDS);
		//cplex.setParam(IloCplex.DoubleParam.TiLim, 60.0);

//		cplex.setParam(IloCplex::Probe, -1);
//		cplex.setParam(IloCplex::MIPEmphasis, IloCplex::BestBound);
//		cplex.setParam(IloCplex::BndStrenInd, 0);
//		cplex.setParam(IloCplex::RepeatPresolve, 0);
//		cplex.setParam(IloCplex::HeurFreq, -1);
//		cplex.setParam(IloCplex::FPHeur, -1);
//		cplex.setParam(IloCplex::NodeSel, IloCplex::DFS);

		cplex.setParam(IloCplex::Threads, 1);
//		cplex.setParam(IloCplex::MIPDisplay, 4);

		// Optimize the problem and obtain solution.
		int status = cplex.solve();
		if (!status)
		{
			env.error() << "No se pudo resolver el modelo. " << status << endl;
			throw(-1);
		}

		env.out() << "Solution status = " << cplex.getStatus() << endl;
		env.out() << "Solution value  = " << cplex.getObjValue() << endl;

		env.out() << "Número de vertices = " << _grafo.getN() << endl;
		env.out() << "Densidad  = " << _grafo.getDensidad() << endl;
		env.out() << "Clique in initial model  = " << branching->getInitialModel() << endl;
		env.out() << "Clique separated  = " << separator->getCliques() << endl;
		env.out() << "Rango separated  = " << separator->getRangos() << endl;
		env.out() << "Rango con pesos separated  = " << separator->getPRangos() << endl;
		env.out() << "Número de nudos B&B  = " << cplex.getNnodes()+1 << endl;
		env.out() << "Número de llamadas a PL  = " << separator->getNLlamadas() << endl;
//		IloNumArray solucion(env, n);
//		cplex.getValues(x, solucion);
//		env.out() << "Values        = " << solucion << endl;
//		IloNumArray vals(env);
//		cplex.getSlacks(vals, constraints);
//		env.out() << "Slacks        = " << vals << endl;
//		cplex.getDuals(vals, constraints);
//		env.out() << "Duals         = " << vals << endl;
//		cplex.getReducedCosts(vals, x);
//		env.out() << "Reduced Costs = " << vals << endl;

		// Genera la soluci�n
//		ret = new DoubleValueArray(n);
//		for (int i = 0; i < n; i++)
//			ret->setValue(i, solucion[i]);

		delete separator;
		delete &x;
		delete heuristicas;
		delete branching;
		delete intera;
	} catch (IloException& e)
	{
		cerr << "Concert exception caught: " << e << endl;
	} catch (...)
	{
		cerr << "Unknown exception caught" << endl;
	}

	env.end();

	return ret;
}

DoubleValueArray* CplexStabSolver::solveLinearRelaxation()
{
	DoubleValueArray * ret = 0;
	IloEnv env;
	try
	{
		IloModel model(env);

		IloNumVarArray var(env);
		IloRangeArray con(env);

		int n = _grafo.getN();

		IloEnv env = model.getEnv();

		//variable definition
		IloNumVarArray x(env);
		for (int i = 0; i < n; ++i)
			x.add(IloNumVar(env, 0, 1));

		//objective function
		IloNumExpr fobj(env);

		for (int i = 0; i < n; ++i)
			fobj = fobj + x[i];

		model.add(IloMaximize(env, fobj));

		//constraints
		IloRangeArray constraints(env);
		for (int i = 0; i < n; ++i)
			for (VertexSetIterator j = _grafo.vecinos(i).begin(); *j < i && j != _grafo.vecinos(i).end(); ++j)
				constraints.add(x[*j] + x[i] <= 1);
		model.add(constraints);

		IloCplex cplex(model);

		// Optimize the problem and obtain solution.
		if (!cplex.solve())
		{
			env.error() << "No se pudo resolver el modelo." << endl;
			throw(-1);
		}

		env.out() << "Valor relajación lineal  = " << cplex.getObjValue()	<< endl;

		for (int i = 0; i < n; i++)
			env.out() << "x[" << i << "]=" << cplex.getValue(x[i]) << endl;

	} catch (IloException& e)
	{
		cerr << "Concert exception caught: " << e << endl;
	} catch (...)
	{
		cerr << "Unknown exception caught" << endl;
	}

	env.end();

	return ret;
}
